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We study the thermalization of the classical Klein-Gordon equation under a u 4 interaction. We 
numerically show that even in the presence of strong nonlinearities, the local thermodynamic equi- 
librium state exhibits a weakly nonlinear behavior in a renormalized wave basis. The renormalized 
basis is defined locally in time by a linear transformation and the requirement of vanishing wave- 
wave correlations. We show that the renormalized waves oscillate around one frequency, and that 
the frequency dispersion relation undergoes a nonlinear shift proportional to the mean square field. 
In addition, the renormalized waves exhibit a Planck-like spectrum. Namely, there is equipartition 
of energy in the low-frequency modes described by a Boltzmann distribution, followed by a linear 
exponential decay in the high-frequency modes. 



I. INTRODUCTION 

The Klein-Gordon (KG) equation describes a wide va- 
riety of phenomena, including both classical wave sys- 
tems, such as the displacement of a string attached to 
an elastic bed [l| , or semiclassical and quantum systems 
based on scalar field theories [2]. In addition to the lin- 
ear dispersive terms, KG models often include nonlin- 
ear terms such as a quartic (it 4 ) potential. For instance, 
modern applications of the classical u KG wave equation 
arise both in models of early cosmology and in ultrarela- 
tivistic heavy ion collisions Q . In the context of a string 
on an elastic bed, the addition of a u 4 potential alters 
the elastic bed from a collection of linear springs to a 
collection of nonlinear ones. 

In applications such as the early universe [§], the long 
time-statistical wave behavior governs the equilibrium 
physics. The thermalization, or equivalently the distribu- 
tion of wave energy throughout the Fourier modes, pro- 
vides a particularly useful description of the equilibrium 
state. As a result, recent thermalization studies for the 
KG equation have been done in both quantum [1, [||, 
and classical 0,0,0] field theories. The studies indicate 
that generic initial wave fields tend to thermalize into 
a state with large quantities of energy in a wide range 
of Fourier modes. Some care however is required when 
interpreting the results of any numerical experiment, as 
there are differences between a finite lattice and a con- 
tinuous differential equation. In this article, we focus on 
characterizing the thermal state of the classical Klein- 
Gordon (it 4 ) wave equation, with particular emphasis on 
managing the ultraviolet s pec trum. 

In previous studies [!, [rjj-llOl] , the thermalization of the 
classical KG equation has been examined with an em- 
phasis on the applications to quantum field theory or 
quantum systems. For example, Boyanovsky et al. Q 
examine the approach to thermalization for out of equi- 
librium initial conditions. Meanwhile, Aarts et al. |6[ 



consider both the thermalization of single fields as well 
as the statistical average of many initial fields (canoni- 
cal ensemble average). Moreover, they also consider the 
interesting case of a strong nonlinearity as a nonpertur- 
bative system. One general trend in the previous work 
is the existence of a local thermodynamic equilibrium 
(LTE) P. 

In our work we study the LTE solutions for a wide vari- 
ety of initial wave configurations. Hence, we evolve each 
solution from a fixed, out of equilibrium initial condition 
and do not consider ensemble averages. Instead, all ap- 
propriate thermal quantities are obtained by averaging 
over time. For instance 



t+AT/2 



A(t')dt', 



(1) 



t-AT/2 
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where AT is the averaging interval. The quantity (A), 
however, is not fixed for all times, but exhibits character- 
istics of an LTE and (for a fixed averaging interval) may 
drift in time. 

In contrast to previous work, which use an effective 
Green's function (two point correlation function) Q, we 
introduce a renormalized wave basis defined by the re- 
quirement of vanishing wave correlations. Such a method 
is described by Gershgorin et al. fl2l. ]l3j. however, their 
work applies to the case of a finite lattice and not a dif- 
ferential equation. A similar recent study [lij ]. carried 
out in the Majda-McLaughlin-Tabak wave system, also 
examines the resulting renormalized dispersion relation, 
and demonstrates how the new dispersion relation can 
effect the dynamics of the wave resonance structure. In 
analogy with previous work, we show that the renormal- 
ized basis exhibits features of a weakly nonlinear system, 
even in the presence of strong nonlinearities. As a re- 
sult, we obtain a simple form for the renormalized wave 
dispersion relation. Here we find that the renormalized 
dispersion relation is not related to the linear dispersion 
relation by an amplification factor, as found in the case 
of a lattice, but rather undergoes an effective mass shift. 
Moreover, we find that the effective mass shift is differ- 
ent than that suggested by the simple (Hartree) nonper- 
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turbative approach. In the classical case of a string on 
an elastic bed, the characterization of the LTE may be 
posed as understanding the long time-averaged sound of 
the string. For instance, what frequencies does the string 
make (dispersion relation), and how loud do they sound 
(energy distribution)? 

This article is organized as follows. We first introduce 
the Klein-Gordon equation and the coordinate transfor- 
mation to renormalized waves. Secondly, we outline our 
numerical experiments, and provide verification that the 
renormalized basis exhibits characteristics of a weakly 
nonlinear system. In addition, we numerically show the 
nonlinear dispersion relation undergoes a frequency shift 
proportional to the time-averaged mean field, S(t): 



1 f 2v 



u 2 (x, t) dx. 



(2) 



Lastly, we study the LTE spectral energy distribution, as 
well as fluctuations about the equilibrium state. 



II. THE KLEIN-GORDON EQUATION 

We study the classical it 4 Klein-Gordon (KG) field in 
1+1 dimensions. We fix the domain [0, 2ir], so the Hamil- 
tonian is 



H = 



2 IT 



l(ul + ul + u 2 ) + ^dx, (3) 



where A is the coupling constant. Letting p = u t so that 
H = H[p, u] is a functional of the fields p and u, one 
obtains the KG equation via, pt 



ML- 

Sp 



Utt 



(4) 



In addition to the conservation of total energy H, pe- 
riodic boundary conditions u(0,t) = u{2it, t) result in a 
second integral of motion, the total momentum 



P 



271 



u x Ut dx. 



(5) 



\u 2 



ju 4 with A > 



We choose the potential V(u) 

0, intentionally as a convex bowl to ensure the absence 
of coherent wave structures or the localized trapping of 
wave energy [ll|. 



energy in each mode. Specifically, the Fourier series for 
p{x,t) and u(x, t) are 



p[x,t) = 
u{x, t) = 



1 



Pk(t)e~ 



u k (t)e 



-ikx 



with summation k over all integer values. Taking p(x, t) 
and u(x, t) as real fields restricts u k = u*_ k and pk = p*_ k 
so that H becomes 



E h 



; (\ Pk \ 2 +u>i\u k \ 2 ). 



(6) 
(7) 



Here E k is the linear energy in mode k, while io\ = 1 + fc 2 
is the linear dispersion relation. In addition to construct- 
ing Fourier mode solutions, one may also make a sec- 
ondary transformation to study the kinetics of interact- 
ing waves 



CLk 



p k - lUkUk 



In this new basis, the equations of motion become 



ia k 



d 

— H(a k ,a* k ), 



(8) 



(9) 



where the Hamiltonian, H(a kl a* k ) also acquires a conve- 
nient form: 



H = y~^fc|afc| 2 - 

k 



(10) 



With the dynamic equations written in the variables 
a k , the linear oscillator solutions take the form a k (t) = 
A k e~ ZUJkt , where the amplitude and phase are contained 
in the complex variable A k . As a result, each wave 
number oscillates with a negative frequency so that the 
Fourier transform of a linear wave is a Dirac 8 function 



2fc(w) 



1 



V2tt J- 
A k 8(w k + lu) 



a k (t)e~ 



At 



(11) 
(12) 



Secondly, the infinite time correlation of any two linear 
waves, including the case of k — ±i, vanish: 



A. The linear case 

Although we study the wave system (0} for strongly 
interacting waves, to motivate the introduction of renor- 
malized waves, we first address @ in the linear case, 
followed by the case of weakly interacting waves. In the 
absence of nonlinearity, A = 0, (0} reduces to the linear 
KG equation. Hence, wave solutions decouple into lin- 
early independent Fourier modes, with conservation of 



\a k ai) 



0. 



(13) 



Hence, linear solutions are uncorrelated waves, which os- 
cillate at a single frequency. 



B. Weakly nonlinear renormalized waves 

In the presence of nonlinearity, A > 0, the waves a k no 
longer decouple into linear oscillators and therefore no 
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longer conserve the energy E k . For instance, the wave 
amplitudes |dfc(£)| do not remain constant but fluctuate 
from the nonlinear interaction. As a result, we seek to 
characterize the nonlinear effect on both the amplitude 
and phase behavior of the waves a k . Before proceeding 
to the case of strongly interacting waves, we first exam- 
ine the case of a weak nonlinearity. By weak nonlinearity 
we mean the maximum amplitude of u is order 1, while 
the coupling strength remains small A <C 1. Moreover, 
with the onset of nonlinearity, especially in the presence 
of strong interactions, the waves a k no longer exhibit the 
linear properties (|12p and (TT3"1) . We therefore follow the 
ideas of Gershgorin et al. [12j and introduce renormalized 
waves Cfe in lieu of (J5J . Renormalized waves form a useful 
basis since they exhibit characteristics of linear waves. 
Namely, renormalized waves have vanishing correlators 
and approximately oscillate with a single frequency. In 
analogy with the linear waves ([5]), renormalized waves 
are defined through a linear combination of the Fourier 
transforms u k and pk, however, linear frequencies are al- 
tered to renormalized frequencies: 



C k = —==(p k - ILUkUk) 



(14) 



The transformation to renormalized waves is canonical 
(up to a factor of i) provided the yet to be determined 
frequencies satisfy uj k > and Cj k = 

We now examine the evolution of renormalized waves 
in the presence of a weak nonlinearity A 1 through a 
perturbation expansion. Specifically, we show that the 
standard procedure of choosing renormalized frequencies 
to eliminate resonant terms yields the same result as 
choosing frequencies to enforce a vanishing correlator for 
waves at ±fc wave numbers, i.e., (c k c_ k ). Written in 
terms of renormalized waves, the dynamic equations for 
Cfc are 



ICk 



ic k 



dH(c k ,c* k ) 



Tklr, 



r 

l,m,n 

- 12c lC * m c* n 8i_ k - 
1 



=fc + (1 - 4) C -fc 
w fc 

^ C * C m C n^k+l+m+n 
— m — n 4QC m C n (5fc_/_ 



(15) 
(16) 



327r(d)fc(Dzdv n ,d;„ 



where Tki mn is the series coefficient, and S x = 1 when 
x = and S x = otherwise. The summation is over 
all integers l,m,n. As is commonly the case in kinetic 
theories, wave behavior is dominated by interacting reso- 
nant terms, for instance, those which force the order-one 
linear oscillators c k on resonance. The nonlinear term 
in (|15l) admits only one such resonance, which comes 
from the term c„c m c ; * = \c n \ 2 Ck when k = n,m = I 



or fc = m, n = I. The presence of only one resonant term 
follows from the fact that the linear dispersion relation 
uj\ = 1 + fc 2 is concave up [16]. For instance, the har- 
monics generated by the products c*c^c*, qc^c*, and 
Qc m c„ cannot oscillate with frequency uj k provided Cj k is 
a concave function of fc. 

In addition to each resonant term in the series c n c m c ; *, 
there is an equivalent term in c/c* n c* containing the vari- 
able c*_ k when n = — fc, I = m or m = — fc, I — n. We now 
remove both the resonant terms and their equivalent c*_ k 
terms from the nonlinear series, and absorb them into 
the coefficients of c k and c*_ k . We do this so that upon 
seeking a small amplitude (A -c 1) solution, we may si- 
multaneous handle both the zeroth, O(l), and first order, 
0(A), corrections to the renormalized frequencies. The 
dynamic equations then become: 



«c fc = g k (t)c k + h k (t)c*_ k + Xr k (t), 



(17) 



h k (t) 



jUfe(i) — 24 ^T mm fcfc | c m | 2 — l2T kkkk \c k \ 



r k{t) — ^ Tkimn 
l,m,n 

- 12 cic*cZ5, 



■m~n"l-k-rn—n 
\2<C n C m Ci Sl+lz — ni—r, 



^cic m c n 5 k —i— rr i—f 



Here the prime in r k restricts the summation to exclude 
all resonant terms having either c k or their equivalent 
terms containing c*_ k . Meanwhile, the time dependent 
coefficients g(t) and h(t) are related by 



9k{t) + h k (t) = Q k . 



(18) 



Thus far, Eq. is exact. 

We now extract approximate solutions to ([T7]) . in the 
limit A 1, for the initial value problem: 



Cfc(O) - C fc . 



(19) 



To remain consistent in the small amplitude approxima- 
tion, we assume the initial field is small with finite energy. 
Hence, each C k < 1 while C k — > as fc — > oo. To obtain 
a solution, we seek an asymptotic series for c k in powers 
of A with a single frequency leading term: 



Cfc(t) = C k e-^ kt + Aci 1} (t) + \'c{:>{t) + 



(20) 



In general, even at zeroth order in A, an arbitrary choice 
of cjfc will couple solutions c± k together through the co- 
efficients g(t) and hit). As a result, the proposed ansatz 
([20]) requires the following two consistency conditions on 
Wfc: 



9k 
h k 



Wfe, 

0. 



(21) 
(22) 
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The first consistency condition (|2Tj) comes from the fact 
that we have chosen c k to oscillate with frequency u k . 
Meanwhile, the second condition (f22|) follows from the 
fact that our ansatz has no dependence on the initial 
value C_fe, or equivalently that Ck decouples from c_fc. 
Relation (fl~8|) however, implies the equivalence of the two 
conditions ([2~Tj) . (j2"2")l . Therefore, choosing uj through (J^TJ) 
as the Ck oscillator single frequency automatically guar- 
antees (|22|) . thereby decoupling the fields Ck and c*_ k . 
Using Eq. (|2"TT) to extract the expression for oj k yields 



UJk 



(0) 



2 V 1 + ^ 



0(A 2 ), 



24 \ ^mmfcfc|C'm|* 



12TfefefejJCfe] . 



(23) 
(24) 



To solve for Cb k to 0(A 2 ), we assume the amplitude in a 
single mode is much smaller than the total averaged wave 
field. Therefore, neglecting the term T kkkk \C k \ 2 results 
in 



u>l = loI + 3\{S) + 0(A 2 ) 
= l + 3X(S)+k 2 +0(X 2 ). 



(25) 
(26) 



Here we have made use of the fact that /i? is proportional 
to the time-averaged value of the wave field 



{S{t)) = {^jJ u 2 {x,t)dx) 



^4^r (|c, - m|2 + |Cm|2) - 



(27) 
(28) 



Note that within the framework of the small amplitude 
ansatz, the approximation (|25[) only remains valid over 
time scales 0(A~ 2 ). As a result, to remain consistent, 
the time average taken in (|2"T|) should be made over, at 
most, comparable time scales, <3(A~ 2 ). 

With the values of ui k solved to 0(A), we may substi- 
tute cf^ = C k e~ lU3kt into the nonlinear term r(t) and 
solve for the first-order contribution to c k (t). Conse- 
quently, c k (t) satisfies the equation 



ui k c k 1] 



(t), 



4 1} (0)=0. 



(29) 
(30) 



Here r[° (i) is the nonlinear term r(t) evaluated using 
C k e~ tUkt . Equation (l29l) therefore describes a 
first-order linear oscillator forced off resonance by r k °^ (t) . 



,(o) 



The solution for c k (t) is 



- klmn 



l,m,n 

12QC*C*e* 



J i 

Cb k + Cui + 



5k+i+ 



m+n 



Wfe — u>i + u> m 4 




4G l C m C n e-« s *+ a > 


- — s 


d k — &i — u m — 






„-w„)t 



n-k- 



k—l — m—n 



LU k + LUl 

De~ lQkt . 



Si+k- 



(31) 



In the solution for c k X \ the prime in the summation 



indicates that there are no terms oscillating with fre- 
quencies ±o)fc. This follows from the fact that we have 
removed the terms in r(t) oscillating at frequencies Oo k . 
Meanwhile, the constant D, from the homogeneous term, 
is chosen to satisfy the initial value (0) = 0. With the 
asymptotic solution (f2"U|) solved to 0(A 2 ), we may cal- 
culate the time averaged correlator, (c k c^ k ), for waves 
±k: 

-i2Cb k t\ | \r<. /„— tcDj.t„(l) \ 



(ckc.k) = CkC^kie- 1 ^) + \C k {e 



(32) 

+ XC-kie-^'ci:') +0(A a ). (33) 

Using formula (JT]) with an averaging time of 0(A 2 ) will 
0(A 2 ). Meanwhile, the long time 



result in 
average (< 



-4lifc* c ( 1 )\ 



picks out the component of c_L with 



frequency oj k . Since c+ k (t) contains no term oscillating at 
frequency uj k , taking a time average over length 0(A~ 2 ) 
will also tend this term to zero: 



(c k c- k ) = 0(A 2 ). 



(34) 



Therefore, to first order in A, choosing the renormalized 
frequencies u k to remove the resonant terms in (|15p is 
equivalent to choosing them to enforce a vanishing cor- 
relator. 



C. Strongly nonlinear waves 

In studying LTE solutions, we are specifically inter- 
ested in the case of strongly interacting waves, A > 1. 
Even in the presence of a strong nonlinearity, one may 
introduce renormalized waves defined by (I14[) . Extend- 
ing the properties of linear waves, the new renormalized 
frequencies Co are then found by imposing that waves Ck 
and C-k remain uncoupled and their correlations vanish 
El El: 



(CkC-k) 



0. 



(35) 



If one has a KG solution u{x, t) for some long time inter- 
val, substituting the definition of Ck into (|35j) yields an 
expression for uj k : 

(36) 



(\Pk 



(Nl 
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Equation (J36J) provides a useful formula for numerically 
obtaining the renormalized frequencies. Despite the fact 
that the renormalized transformation defined by Eqs. 
(H3J) and ([3"o]) is always mathematically possible, there is 
no guarantee that the waves Ck stand out as a useful co- 
ordinate system. The waves Ck do however form a useful 
basis provided the Fourier modes exhibit narrow band 
single frequency oscillations. When Uk does oscillate with 
roughly one frequency, the formula (|36[) identifies uj^ with 
the frequency of oscillation. In the case of the KG equa- 
tion with a strong u A interaction, we do in fact find that 
the waves Ck form narrow band oscillators with a central 
frequency We verify these wave properties in our nu- 
merical experiments, implying renormalized waves are a 
natural basis to study the LTE spectra. We demonstrate 
in our numerical experiments that in the limit A ^> 1, 
the renormalized frequencies are very well approximated 
by the dispersion relation 



,~.2 



1 + 2.59A(S*) + fc 2 



(37) 



Here we refer to 2.59A(5) 



u)\ as the mass shift. 



In the case of large A, the nonlinear time scale may 
be estimated by comparing the size of terms in the KG 
equation. Specifically, assuming the field u(x, t) is or- 
der one, the time derivative will balance the nonlinear 
term, u tt ~ Au 3 , provided that the characteristic non- 
linear timescale tq ~ A~ 1//2 . Note that this estimate 
only holds for large A, or when A dominates the linear 
dispersive term in the KG equation. In our numeric ex- 
periments, we find that LTE solutions exhibit 0(A -1 / 2 ) 
as the natural fundamental period, while earlier work fo- 
cusing on the analytic development of short time KG so- 
lutions [l7| further verifies A -1 / 2 as the strong amplitude 
time scale. Consequently, we define LTE time-averaged 
quantities (A) over many periods of the time scale tq: 



t+AT/2 

(A)(t)= A Um o ^ / A(t')dt'. (38) 



III. THE LOCAL THERMODYNAMIC 
EQUILIBRIUM 

A. Equilibrium properties 

In our numeric experiments, we study LTE solutions 
for a wide range of initial conditions. Specifically, we fo- 
cus on initial conditions set over a low band of Fourier 
modes, and vary both the number and amplitudes of 
modes excited. In analogy with experiments on the 
Fermi-Pasta-Ulam (FPU) model [H-|2l|, we find that 
there is a generic transient stage where energy in the 
initial modes is redistributed to higher ones. For large 
nonlincarities, A > 1, the transient stage appears to oc- 
cur on a timescale at least as fast as r ~ A" 1 / 2 . Upon 
redistribution, the energy spectrum settles into an LTE, 
where a majority of the energy is shared within a finite 
band of lower modes. Here the LTE persists for time 
scales much longer than tq and appears similar in nature 
to the intermediate metastable state realized by similar 
experiments [l9l - |2~i| in the FPU lattice. 

In a recent work concerning the FPU lattice, [2l[ the 
authors show that initial conditions, specifically those 
with Fourier components having random phases versus 
coherent phases, can alter the long time behavior of the 
solutions. Although most of our test cases under consid- 
eration have initial conditions with coherent phases, we 
also ran experiments with uniformly distributed random 
phases. We found that the effect of the phases did not 
effect the onset, or characteristics of the LTE. 

Since the characteristics of a local thermodynamic 
equilibrium appear to be robust for the various initial 
conditions that we test, we first focus on describing in de- 
tail the generic long time behavior of one solution u(x, t) 
to Eq. @ , and defer a description of the general trends 
for the following subsection. 

Rescaling the field u(x,t) — > \~2 U (x,t) effectively sets 



the coupling constant to unity: u u 



-u—u A . Hence, 



without loss of generality, we take A = 1 and control the 
coupling strength through the initial Fourier amplitudes. 
For the test case under consideration, we naturally ini- 
tialize u(x, 0) via Fourier components [22[ to 



t-AT/2 



Uk = 0, 



(39) 



Our numeric experiments try to capture solution prop- 
erties at large wave numbers, especially within the expo- 
nential decay of the spectrum. We therefore use a pseu- 
dospectral method to maintain high numeric accuracy in 
both space and time. The pseudospectral method con- 
sists of exact spectral propagation for linear terms, cou- 
pled with a Richardson extrapolation [l8[ algorithm for 
the propagation of the nonlinear terms. For smooth solu- 
tions, it, the spectral basis ensures high spatial accuracy, 
while the Richardson extrapolation guarantees high ac- 
curacy in time. Appendix B contains a more detailed 
discussion of the numerical algorithm. 



_ J 2.08 < \k\ < 50 
Uk ~\0 k = 0, \k\>50 

Although each Uk is O(l), the initial mean squared field 
amplitude is in fact strong: S = 68.7. The total energy 
and momentum are H ~ 6.7 x 10 5 and P = 0. 

We obtain the renormalized waves and study the KG 
thcrmalization in two steps. First, we evolve the field 
u(x,t) for a long period of time, i.e., from T = to T\. 
We choose T\ long enough so that the field amplitudes 
complete the initial transient stage and redistribute 
their energy into higher modes. More specifically, follow- 
ing [19h21| one may introduce M(t) as the number of 
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FIG. 1. Spatiotemporal spectrum of renormalized waves 
taken over the interval (153,230) for initial data (|39|l . The 
shaded plot shows | Cfc (oj) | 2 while solid line corresponds to the 
dispersion relation defined by 




100 



cal- 



FIG. 2. Power spectrum of renormalized waves |ci;(u;)| 2 
culated over the interval Ti = (153, 230) for initial data 
The waves stay localized around — u>k while the oscillations 
around ujk have been completely removed. From right to left, 
the peaks correspond to k = 5, 20, 40, 60. 



thermodynamic equilibrium and both averages and fre- 
quencies Qk do depend on the interval (T 1 ,T 2 ). Despite 
the fact that the wave equation only achieves a LTE, the 
LTE does exhibit generic properties. 

To verify that renormalized waves oscillate with ap- 
proximately one frequency, and that the frequency is in 
fact (life, we calculate the spatiotemporal transform of Ck 
over interval times (Ti,T2). Figure (JXJ) shows the spa- 
tiotemporal transform |c\.(w)| 2 for an early time interval 
Ti = 153, T 2 = 230. Here the choice of times T x and T 2 is 
somewhat arbitrary. Indeed we also obtain similar plots 
for other time intervals, such as T\ = 767, T 2 — 843 or 
even Ti = 153 and T 2 = 843. To calculate the spatiotem- 
poral transform, we first compute il>k using (|36[) averaged 
over (Ti,T2) to obtain Cfc(i). For each k, we compute 
the power spectrum |cfc(w)| 2 . We find that the energy 
in each Ck(u>) is sharpely peaked, with a small amount of 
background noise. The background noise can be removed 
by treating each Cfc(i) as a stochastic signal, and averag- 
ing multiple power spectra. For instance, to remove the 
background noise, we divide the interval (Ti,T 2 ), which 
contains 25000 data points, into eight partitions. We 
calculate the power spectrum of each partition indepen- 
dently, and as shown in figure average the results 
together. In addition to the spatiotemporal transform, 
the thick line in figure ([1]) is the dispersion relation tuk 
obtained via (|36[) . The dispersion relation ujk corresponds 
very well with the localized magnitude of |cfc(w)| 2 , indi- 
cating the waves Ck effectively oscillate with frequency 
— cDfc. The faint line observed at Cjk is a small residual of 
the coupling found between Ck and c*_ k . Since the plot 
is on a log shading scale, the coupling effect is in fact 
several orders of magnitude smaller than the amplitude 
found at — Q).. The removal of the coupling found be- 
tween Ck and c*_ k is further verified by the absence of a 
peak in the power spectrum at ujk seen in figure 

In addition to numerically computing ujk, we find that, 
in the case of strong nonlinearities, the frequencies always 
satisfy a shifted Klein-Gordon dispersion relation of the 
form 



uj% 1 + 2.59A(S) + r 



(40) 



modes which share most (i.e., 90%) of the energy. Dur- 
ing the transient stage, the value of M(t) rapidly, over a 
time scale at least as fast as A -1 / 2 , increases to a plateau. 
We use the onset of a plateau region in M(t), and hence 
the stable sharing of energy between a finite band of low 
wave number modes, as an indication for an LTE. 

Once the field has thermalized, we then further evolve 
the modes Uk and pt over some time interval (Ti,T 2 ) to 
obtain a description of the thermodynamic equilibrium 
state. Specifically, we recast the fields Uk and pk into 
renormalized waves by obtaining the appropriate time 
averaged quantities, ( | "tife 1 2 ) and (|Pfc| 2 ), followed by 
and Cfc using Eqs. (|14[) and (|36p . In agreement with 
previous studies, we find that waves achieve only a local 



Again S is the mean squared averaged field computed 
over the same time interval (Ti,T 2 ) as <Dfe. For exam- 
ple, figure pip shows the dispersion relation for several 
different initial conditions. Moreover, as a result of the 
LTE, both the value of S and therefore Qk drift over long 
times. Figure © verifies the relation (T4"0"|) on two sep- 
arate time intervals, with (S) = 67.4 on (153,230) and 
(5) = 60.24 on (767,843). Figure g]) also shows a plot 
of S(t) highlighting the two averaging intervals. The de- 
viation in (5*) between the two intervals corresponds to 
a drift of roughly 5% in the smallest frequency (Do 

Since the spatiotemporal transform verifies that the 
waves Cfc form narrow band oscillators centered at fre- 
quency Uk , we may introduce the effective energy [l3j in 
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FIG. 3. Frequency drift of renormalized waves. As the aver- 
age S drifts in time, so do the renormalized frequencies. The 
top and bottom curves represent the spectrum of renormal- 
ized waves over the time intervals (153,230) and (767,843), 
respectively, for initial data p9[) . The dotted curves show the 
frequencies obtained by equation (136[) . while the solid lines 
represent fits obtained using (|40|l . The quantity S is averaged 
over the same interval used to calculate the wave spectra. 
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FIG. 4. Plot of the mean field S(t) for initial data ([39]). The 
dark curve represents a local time average (S) [over the in- 
terval (t — 14, t + 14)], and demonstrates a drift over time. 
The double arrows highlight the intervals h = (153, 230) and 
h = (767, 843). 

Fourier mode k in analogy with a linear oscillator 

(E k )=±({\p k \ 2 )+u>l{\u k \ 2 )) (41) 

= <M 2 )- (42) 

Here the definition of the renormalized frequency over a 
given interval implies the waves on average satisfy the 
virial theorem, equally splitting the kinetic and potential 
energy. The quantity (|j>fc| 2 ) therefore is a measure of the 
approximate energy in mode k, independent of renor- 
malized frequency. Although the energy E k describes 
the modal distribution of energy, one should note that 
~Y^u{E k ) is not a conserved quantity and that the values 
of (Ek) exist only in an LTE state, and weakly depend on 



the interval (Tx,T2). The energy E k can also be related 
to the renormalized wave amplitudes via 

£* = f (<k| 2 } + <|c- fc | 2 )). (43) 

Consequently, in the special case when Uk = u-k, such 
as the initial data (13"9"|) . then Ek = ojk(\ck\ 2 ) can also be 
used as a measure of energy in mode k. 

In general, over the thermalization stage (0,21), so- 
lutions initialized to lower Fourier modes leak out into 
the higher modes. When viewed on a log scale, the en- 
ergy spectrum evolves into a flat distribution in the low 
Fourier modes accompanied by a linear exponential de- 
cay in the high modes. For example, figure ([5]) shows 
the energy spectrum (Ek) averaged over the long time 
interval (152, 843). As solutions propagate through time, 
the spectrum retains a straight exponential decay, how- 
ever the slope of decay may drift mildly (e.g., 10%) over 
time. We also calculate the sum J2k(Ek) ~ 6.8 x 10 D 
which differs from the exact energy H ~ 6.7 x 10 5 by 
« 2.1%. Qualitatively, the spectrum appears very similar 
to the one predicted by the Planck blackbody distribu- 
tion: |cfc| 2 cx (e l3uik — 1) _1 . The presence of a Planck- like 
spectrum arising from a classical system, such as a sys- 
tem of weakly coupled oscillators, is also discussed by 
Carati and Galgani (23|. For our data, the fit is only 
heuristic in that we can not simultaneous match the flat, 
long wave spectrum as well as the short wave exponential 
decay with a fixed temperature. 

In addition to studying initial conditions (I3"9")) for TV = 
2048 modes, we also performed several convergence tests 
with varying grid sizes TV = 512, 1024, and 4096 and dif- 
ferent time steps At. We also found that numeric solu- 
tions do not depend on grid spacing provided the power 
spectrum remains exponentially suppressed over the in- 
tegration times. For instance, TV — 4096 with initial 
(l3"9l yields identical solutions to TV = 2048. Meanwhile, 
TV = 512 develops an energy cascade into the ultravio- 
let spectrum, accompanied by the onset of equipartition 
of energy. Hence, for the initial conditions (|3"9"|) . taking 
TV = 2048 yields a consistent solution to the PDE over 
the integration times < T < 843, while inconsistencies 
develop for smaller TV < 1024. 



B. Kinetics and fluctuations of the LTE 

In addition to the power spectrum and renormalized 
dispersion relation, we also study the kinetic behavior of 
the LTE as well as fluctuations about equilibrium values. 
Since the frequency shift tracks the spatial average 5, we 
plot S as a function of time to gain an understanding 
of the dispersion relation evolution. Figure (@| shows 
the value S over our integration time, along with a local 
time average value of S. We take the local time average 
S over intervals long enough for the waves c k to exhibit 
numerous oscillations. Although the waves Ck appear well 
defined for any fixed time interval (7\, T2), there appears 
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no well denned equilibrium, or more precisely if even the 
limit liixiT-yoo y Jq S dt exists. As shown in figure (U]), S 
generally appears to decay, however drifts in time. 

To further characterize the nature of the LTE, we ex- 
amine fluctuations by plotting probability distributions 
of the squared field strength S and energy spectrum for 
different time intervals. We find that for a sufficiently 
long time interval, S acquires a Gaussian probability dis- 
tribution with varying mean and width. Figure ([S]) shows 
the probability distribution S, as well as a Gaussian fit, 
over the interval (153,230). Although not shown, the 
interval (767, 843) also exhibits a Gaussian probability 
distribution with a different mean. Only at very short 
time intervals, does S not widen out to a Gaussian. 





FIG. 5. Energy spectrum of renormalized waves (on log scale), 
Ek = (Jpfc| 2 ). Spectrum is qualitatively similar to the distri- 
bution derived by Planck, equipartition of energy in the low 
modes, exponential decay in the high modes. Data is averaged 
over the time interval (153,843) for initial data (|39l) . 



We also study the probability distributions (PDF) for 
wave amplitudes or particle numbers \ck\ 2 - The distribu- 
tion of these variables are related to the energy spectrum 
since one can identify cD/c|cfc| 2 as the energy in mode fc. 
Using the initial conditions ([39]) we plot the distribu- 
tion of |cfc| 2 over the interval (767, 843) for different wave 
numbers. We find that low wave numbers, fc < 200, 
which achieve equipartition of energy, exhibit a classi- 
cal Boltzmann distribution. For example figure ([7]) il- 
lustrates the distribution for mode fc = 2 along with an 

a i 1 2 

exponential fit of the form e _p|Cfe . Meanwhile, the expo- 
nential distributions found in wave modes near the nat- 
ural spectral cutoff, fc ss 200, start to shift their peaks 
away from |cfc| 2 = 0, as seen in figure |8]). Lastly, there 
does not appear to be a consistent distribution of energy 
at large wavenumbers fc > 200. Specifically, adjacent 
modes fc and fc + 1 may exhibit very different distribu- 
tions. Despite the lack of a unifying distribution, many 
large wave modes do exhibit a multipeaked distribution. 
Figures ((9]) and (fT0|) show the amplitude distribution for 
fc = 510 and fc = 520, respectively. 



FIG. 6. Histogram of mean field S over time interval 
(153, 230) for initial data ([39]). The mean field S is well fit by 
a Gaussian distribution with mean 67.4. 



IV. GENERAL TRENDS 

In the following section we study numerical solutions 
for a wide range, roughly 40 trials, of initial data. In 
our experiments, we evolve initial data for a fixed time 
T\ = 767, after which we calculate quantities character- 
istic of the LTE. We test three generic sets of initial con- 
ditions shown in Tables (J]), (JTXJ) , and (|IIip. Specifically, 
the tables show data for the average (S), the frequency 
shift Qq obtained by a best fit to the renormalized dis- 
persion relation, the variance of 5 1 , the exponential slope 
in the power spectrum and the classical analog of the 
particle number. Here the slope of the power spectrum 
refers to the slope obtained by a least squares linear fit 
to ln(Ek) over the large wave numbers fc. Physically, the 
slope corresponds to one measure of the inverse temper- 
ature /3 = (fcf,T) _1 from a Plank spectrum. In addition, 
the particle number is defined by 



(N)=J2(\ck\ 2 )- 



(44) 



We include the values of (N) as they often play a role in 
quantum mechanics. The renormalized dispersion rela- 
tion fits are also verified in figures (ITT1) . (fl~2"]) and p^|) . 

In our first test trials, we fix the initial shape u(x,0) 
by setting Uk = C a constant over the first 50 modes 
(0 < |fc| < 50). Altering the initial constant C varies the 
nonlinear wave strength. Phenomenologically, trials for 
varying coupling constants exhibit behavior identical to 
the LTE described in the previous section. For instance, 
figure (jlll) shows the renormalized dispersion relation at 
different coupling strengths. 

The second set of tests fix the coupling strength A = 1, 
and varying the initial conditions so that energy H and 
momentum P = remain constant. The primary goal is 
to determine whether drastically different LTE solutions 
arise from initial conditions with the same energy and 
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FIG. 7. PDF for renormalized wave at small wave number, 
(k — 2) | C2 1 2 . Data are taken over the time interval (767, 843) 
for initial data (|39p . The curve represents an exponential, 
Boltzmann distribution. 
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FIG. 8. PDF for renormalized wave for mid ranged wave 
number (k — 209) |c209 1 2 ■ Data are taken over the time in- 
terval (767, 843) for initial data (J39J). Mode k = 209 is on the 
edge of the spectrum between equipartition in the low modes 
and exponential decay in the high modes. The decrease in 
probability at low amplitude shifts the peak to the right. 



momentum P, H. Table (JTTJ) shows LTE trends for initial 
data u k = 0.380,0 < |fc| < 150, u k = 0.648,0 < |fc| < 
100, and u k = 2.08,0 < |fc| < 50. The data show that 
the renormalized frequencies, mean field strength, and 
exponential slope vary dramatically over the initial data. 

Lastly, the third set of initial data fixes A = 1 and 
varies both the initial conditions iik and u k so that the 
total momentum P ^ 0. The goal in this case is to test 
whether moving waves alter thermalization phenomenol- 
ogy. In general, the wave field thermalize into a renor- 
malized wave basis with the characteristic dispersion re- 
lation and Planck-like spectrum. As shown by the 
data, moving waves exhibit trends identical to station- 
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FIG. 9. PDF for renormalized wave for large wave number 
(k — 510) |c5io | 2 . Data are taken over the time interval 
(767, 843) Higher wave numbers exhibit a variety of behavior. 
In the k = 510, wave number, the amplitude jumps around 
between three peaks. 
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FIG. 10. Distribution of renormalized wave for large wave 
number (k = 520) | C520 1 2 ■ Data are taken over the time inter- 
val (767, 843). 



ary ones. 

V. CONCLUSIONS 

We numerically study the long time behavior of the 
classical Klein-Gordon equation with a strong it 4 nonlin- 
ear interaction. By introducing a renormalized wave ba- 
sis, we show the system exhibits characteristics, locally in 
time, similar to a weakly nonlinear system. Specifically, 
the renormalized waves remain uncorrelated and form 
narrow band oscillators centered around one frequency. 
In addition, the renormalized waves, in their LTE state, 
achieve a renormalized dispersion relation described by 
Eq. pi)]) . Here applies to the case of a strong nonlin- 



TABLE I. Different coupling constants. The first column 
shows the initial value for Uk over modes < < 50. Aver- 
ages are taken over the time interval (767, 843) . 



Initial Uk 


H 


P 


(S) 


~2 
UJ 


Var{S) 


-Slope 


(AT) 


0.69 


26607 





8.9 


25.02 


1.9 


0.035 


756 


1.39 


178502 





29.3 


77.7 


4.6 


0.023 


3860 


2.08 


671909 





60.3 


154.0 


9.9 


0.015 


10550 


2.77 


1867200 





89.5 


232.8 


10.8 


0.011 


21175 



TABLE II. Trials with fixed energy H ~ 1.8 X 10 5 and mo- 
mentum P = 0. Initial conditions are taken as constants for 
Uk over the first 50, 100, 150 modes. Averages are taken over 
the time interval (767, 843). 



Num. modes 


Initial Uk 


(S) 


- 2 


Var(S) 


-Slope 


(AT) 


50 


1.39 


29.3 


77.7 


4.6 


0.023 


3860 


100 


0.65 


18.1 


48.0 


2.8 


0.017 


2450 


150 


0.38 


12.1 


33.3 


2.2 


0.012 


1760 



ear coupling strength, but appears qualitatively similar, 
with a different constant, to the one found in the weakly 
nonlinear analysis (1251) . In addition, the mean field (S), 
and subsequently the nonlinear dispersion relation may 
drift as much as 10% over long times, e.g., timescales 
T 3> However, the LTE renormalized dispersion 

relation (T4"0"|) still holds over any time interval (T 1 ,T 2 ). 

We also find that fluctuations about the LTE are de- 
scribed by several characteristic distributions. For strong 
nonlincarities, the probability distribution of the mean 
squared wave field, S, appears as a Gaussian. Meanwhile, 
for low wavenumbers, the amplitude of the renormalized 
waves \ck\ 2 exhibit a Boltzmann distribution, while for 
larger wave numbers, the probability distributions for the 
wave amplitudes no longer form an exponential decay. 
As a result, there is no longer equipartition of energy 
throughout the large wavenumbers. 

In the current study, certain results, such as the exact 
form of the nonlinear dispersion relation, depend explic- 
itly on the u potential. In future work, it would be 
interesting to determine whether similar relations hold 
for generic potentials. 



TABLE III. The trials have data of the form Uk — A withpj, = 
ii k = i-BVl + k 2 , ul = il-k for < \k\ < 50. Trials 1 and 
2 correspond to amplitudes A — B = 0.69, 1.38, respectively. 
Trial 3 corresponds to amplitudes A — 1.38, B — 2A — 2.76. 
Note, we choose pk oc ujk to mimic the initial conditions one 
would require for traveling linear waves. 



Trial 


H 


P 


(S) 


- 2 
LO 


Var(S) 


-Slope 


(AT) 


1 


47208 


41180 


12.2 


33.5 


1.7 


0.021 


1218 


2 


260900 


164700 


33.9 


89.2 


4.6 


0.016 


5047 


3 


508120 


329420 


51.1 


136.5 


7.6 


0.015 


8664 
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FIG. 11. Renormalized dispersion relation for different cou- 
pling strengths over time interval (767,843). Initial con- 
ditions are taken as a constant over the first 50 modes. 
The curves (bottom to top) have initial Fourier amplitudes 
u k = 0.69, 1.39, 2.08, 2.77. The dotted line corresponds to 
the exact renormalized frequencies from (|36p . while the solid 
line is a fit with the numerically calculated dispersion relation 
Co 2 k = 1 + 2.59A(S) +k 2 . 




FIG. 12. Renormalized dispersion relation for initial condi- 
tions with fixed energy H ~ 1.8 x 10 5 and total momentum 
P = 0. Data are taken over the time interval (767,843). 
From bottom to top, the curves correspond to initial data: 
u k = 0.380,0 < \k\ < 150, u k = 0.648,0 < |fc| < 100 and 
u k = 1.38,0 < < 50. 
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FIG. 13. Renormalized dispersion relation for initial condi- 
tions in first 50 modes, with non-zero total momentum. Data 
are taken over the interval (767,843). From bottom to top, 
the curves correspond to integrals of motion: (H, P) ~ {(4.7x 
10 4 , 4.1 x 10 4 ), (2.6 x 10 5 , 1.6 x 10 5 ), (5.1 x 10 5 , 3.3 x 10 5 )}. 



results [25|, |26J outline the absence of a Rayleigh- Jeans 
divergence in classical field theories. Despite the absence 
of a theoretical divergence, the finite truncation in any 
numeric experiment requires a careful interpretation of 
the results. For instance, numerically solving a PDE over 
very long computational times is mathematically identi- 
cal to integrating a discrete lattice with many points. 
More specifically, by truncating the PDE to a finite lat- 
tice with N modes, nonlinear terms may artificially intro- 
duces aliasing effects that remain absent in the PDE. In 
the case of u 4 nonlinearity, upon discretization, Fourier 
modes (k, I, m, n) satisfying the relation k+l = m + n+N 
strongly interact. These aliasing effects, which are ab- 
sent in the PDE system, provide a mechanism for energy 
transport from low to high Fourier modes. As shown by 
De Luca and Lichtenberg [27j, such interactions are re- 
sponsible for the equipartition of energy in a lattice. As 
a result, to distinguish u(x,t) as a solution to the PDE 
and not a discrete lattice approximation, one requires 
that the numerical solutions converge, in an appropriate 
norm, to the continuous field as the number of lattice 
points go to infinity. 
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Appendix B: Numerical Method 



Appendix A: Lattice versus a Partial Differential 
Equation 

In 1965, Fermi, Pasta and Ulam 24] (FPU) performed 
a numerical experiment to show that a weak nonlinear- 
ity would thermalize a discrete lattice. To their surprise, 
they discovered that the lattice did not thermalize to the 
expected equipartition of energy, but in fact the energy 
returned periodically to its initial configuration. In the 
flurry of work following FPU, numeric results exposed 
that the quasiperiodic solutions were in fact dependent 
on initial data and nonlinear coupling strength. For suffi- 
ciently strong nonlinearities, the FPU lattice does achieve 
equipartition of energy [l9j . 

In contrast to the finite dimensional lattice, wave equa- 
tions contain an infinite number of modes. In the context 
of wave thermodynamics, the difficulties with infinitely 
many modes has been known for over 100 years, dat- 
ing back to the Rayleigh- Jeans paradox and ultraviolet 
catastrophe. In particular, the assumption of equiparti- 
tion of energy implies that the total wave energy diverges 
for any fixed thermodynamic temperature. Although this 
difficulty was resolved by Planck for the problem of ra- 
diation, the thermalization of a partial differential equa- 
tion still requires one to work with an infinite number 
of modes. The partial differential equation (PDE) limit, 
which roughly speaking takes the number of discrete lat- 
tice points N — > oo while keeping energy, E and volume 
V constant, is also radically different from the thermody- 
namic limit which takes N — > oo at the same rate as the 
system size V — » oo. One consequence of this difference 
is the absence of equipartition. For instance regularity 



In the following appendix, we describe in detail our 
numerical method for evolving solutions to the Klein- 
Gordon equation. As mentioned in the body of the pa- 
per, we use a pseudo-spectral method, which combines a 
spectral propagation for the linear terms with a Richard- 
son extrapolation for the nonlinear ones. 

To march Eq. Q through time, we convert the PDE 
into an equivalent first-order system: 



U t - iL(U) = N(U) 



(Bl) 



where U = (u,w) T , L and N are linear and nonlinear 
operators respectively. We take coordinates for u and w 
along the characteristics of (U]) so that L and N are 



L 



N = 




(B2) 



Using matrix exponentials, we may integrate the linear 
term in (|B1[) to obtain an expression for the propagation 
over a small time h: 



-iLt 



U(t + h) = e lhL U{t) + e 



N(V), 



t+h 



•(*-«) N[U(s)} ds, 



U{t + h) « e* hL U{t) + hN[V(t)}. (B3) 

When written in Fourier space, the differential operator 
e ihL decouples into 2x2 matrices acting on each Fourier 
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component of U. We therefore use Fourier modes and a 
fast Fourier transform when evaluating the matrix expo- 
nential. Meanwhile, we adapt the Richardson extrapola- 
tion algorithm for propagating ODE's 18], to handle the 
nonlinear term in (|B3I) . 

The Richardson extrapolation routine is a method for 
evaluating numerical limits. In our case, starting with 
U(t) at time t, we seek to obtain the solution XJh(t + At) 
at a fixed time At later, in the limit h — > 0. Once we 
have extrapolated the solution a time step At, we re- 
peat the process M times to obtain a solution U(T) at 
time T — MAt. To accomplish the extrapolation over 
a step At, we truncate our system to N — 2048 spec- 
tral modes and a spatial stencil Ax — 2 ^-. We also fix 
At = 0.5Ax < 1. The limit, \i-m h ^ V(t + At) can now 
be thought of as the simultaneous limit of a finite num- 
ber of variables. Following the Richardson extrapolation 
routine, we integrate solutions over a time At to obtain 
XJh (t + At) using successively smaller values of an inter- 
mediate step h. For instance, taking ho — At in Eq. (|B3[) 
yields a one step Euler approximation Uh (t + At). Tak- 
ing hi = At/2, then requires 2 evaluations of (IB3|) to ob- 
tain U/jj (t+At). In general, the evaluation of U^i+At), 
requires one to divide the time step At into / pieces of 
length h — At jl, followed by integrating (|B3[) I times. 
Provided the nonlinearity in (|B3[) is an analytic function, 
and the solutions u are smooth, a polynomial extrapola- 
tion of the sequence UV , for successively smaller values 
of hj, will converge with error 0(At 2r+1 ). Here r is the 
number of sampled values hj, < j < r. In our numer- 
ics, we have r — 7, which gives us a relative accuracy of 
Au~O(10~ 12 ). 

Over the PDE integration times, we track numerical 
errors by estimating the absolute errors, and recording 
the two integrals of motion H and P. From the Richard- 
son extrapolation routine we obtain an L 1 error estimate 
on the propagated solution for each time step At. Sum- 
ming the L 1 errors at each step then provides an absolute 
bound on the L 1 error between the exact and numeric 
solutions. As a consistency check, we also calculate the 
integrals of motion, and find that they remain constant 
to within roughly one part in 10 12 over each time step 
At. In the worst test cases, the integrals of motion de- 
viate to one part in 10 6 after 3 x 10 5 time iterations. 
Meanwhile, the summed L\ error remains bounded to 
roughly one order of magnitude less than the error in the 
integrals of motion. Lastly, to guarantee the u A nonlin- 
earity does not introduce aliasing effects, we restrict our 
initial data so that Fourier amplitudes at wave numbers 
kmax ~ 600 remain exponentially small over the entire 
integration time. 



Appendix C: Periodic Klein-Gordon Solutions 

In the main body of the paper, we evolve Klein-Gordon 
solutions for initial data with energy in a large number 
of Fourier modes, and show that they thermalize into 



an LTE. Although these LTE solutions appear to arise 
from generic initial data, the Klein-Gordon equation also 
admits traveling wave, exactly periodic solutions [l|. For 
instance, initial data starting on a periodic solution will 
not thermalize, but will remain periodic for all time. In 
this section, we investigate these periodic solutions and 
extract the large amplitude limit A > 1. We show that 
solutions with spatial period 2ir/k oscillate at a frequency 
uj\ w 0.95 + 1.57A(m 2 ) 27 t + k 2 , where (u 2 ) 2 7r is defined 
below. 

Starting with the Klein-Gordon equation 

u u - u xx + V'(u) = 0, (CI) 
V{u)= 1 -v? + ±u\ (C2) 

a periodic, traveling wave solution u, with velocity v, has 
the form 



u{9) = u{6 + 2tt), 
9 = x — vt. 



(C3) 
(C4) 



Substitution then yields an integrable ODE 

V'(u) = 0, (C5) 
(C6) 



{v 2 - l)u" 



~(v 2 -l)(u') 2 + V(u) = V(u m ). 

Here V(u m ) is the constant of integration, where u m 
is the maximum amplitude of the field u. Provided 
v 2 > 1, the Eq. (|C6|) describes a nonlinear oscillator 
for the variable u. Therefore, for a fixed maximum am- 
plitude u m , only specific values of v will yield periodic, 
u(6 + 2tt) = u{9), solutions. Namely, these values of v 
correspond to solutions u which oscillate k times over the 
domain < 9 < 2tt. Hence, v must satisfy 



Ah 



2- 1 ' 2 du 



2n 



y/V{u m ) - V(u) x/^T 



For brevity, introduce the effective mass 
2 r Um 2- 1 / 2 du 



fiumr 1 = 



Wo y/V(u m ) - V(u) 



which yields 



v 2 = l 



f 2 (u m ) 

k 2 ' 



The frequency in time w then becomes 
Lu k = kv, 



u k = k^l + f 2 (u m )/k 2 

= V k2 + P{ u m)- 



(C7) 
(C8) 

(C9) 

(C10) 
(Cll) 
(C12) 



Thus far, the nonlinear frequency relation (IC12I) is gen- 
eral in the sense that one only requires the potential V(u) 
to be monotonic and symmetric, V(—u) = V(u), in u. In 
the original variables x and t, the fcth traveling wave so- 
lution to (|C6p , Uk (9) , is periodic in space x — > x + =£■ and 
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time t t + — . For the special case of a u 4 nonlinear 
potential, the integral f(u m ) and solution u(9) may be 
written down explicitly in terms of elliptic functions [HI ■ 
As a result, we may extract the asymptotic behavior for 
large fields Aw^ 3> 1: 



- 0.71 78 \u 2 m 



1.0458 + 



Au 2 



(C13) 



To compare the frequency shifts for the periodic, trav- 
eling wave solutions to those found in the renormalized 
waves, we recast the maximum field amplitude u m in 



terms of the square averaged field: 

(w 2 )2tt 



1 

~ 2tt\ 


r 2rr 


/ u du 
'o 


k 


1 9 i 

/ u 2 d6» 
'o 


4fc 
~~ 2tt. 


/ ix^ — dii 


4fc / 


'V 2 - 1^-1/2 


~ 2ttV 




Jo 


m 2 d/z 


r 

Jo 


m d/i 



u d^t 



(C14) 
(C15) 
(C16) 
(C17) 
(C18) 



where we have used the ODE (|C6[) to replace -4^ with a 
function of m. Here d[i is the measure defined as 

du 



dfx = 



^V(u m ) - V{u) 



(C19) 



In the large amplitude limit, the averaged field scales 
with the peak field amplitude u m 



(w 2 ) 27r - 0.456947u 2 



0.061347 
A 



O 



1 



X 2 u 2 



(C20) 



Hence, in terms of the averaged squared field (it 2 ^, the 
frequency shift is 



f{u m ) - 1.5708A(u 2 ) 27r + 0.9494 + O 
col w 0.95 + 1.57\(u 2 ) 27r + k 2 . 



1 



Au 2 , 



(C21) 
(C22) 



As with the renormalized wave solutions, the periodic 
solutions Uk(0) exhibit a mass frequency shift propor- 
tional to the averaged field X(u 2 ). In addition, like the 
renormalized waves, the leading term in the large ampli- 
tude shift does not depend on the number of oscillations 
k, however, does admit a different leading coefficient of 
1.57 as opposed to 2.59. 
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